Typeset with jpsj2.cls <ver.l.2.2t» 



Full Paper 



O 



00 



^1- 



Finite- Temperature Transition in the Spin-Dimer Antiferromagnet 

BaCuSi 2 6 

Yoshitomo Kamiya, Naoki Kawashima and Cristian D. Batista 1 

Institute for Solid State Physics, University of Tokyo, Kashiwa, Chiba 227-8581, Japan 
1 Theoretical Division, Los Alamos National Laboratory, Los Alamos, New Mexico 87545, USA 



We consider a classical XY-like Hamiltonian on a body-centered tetragonal lattice, focusing 
on the role of interlayer frustration. A three-dimensional (3D) ordered phase is realized via 
thermal fluctuations, breaking the mirror-image reflection symmetry in addition to the XY 
symmetry. A heuristic field-theoretical model of the transition has a decoupled fixed point 
in the 3D XY universality, and our Monte Carlo simulation suggests that there is such a 
temperature region where long-wavelength fluctuations can be described by this fixed point. 
However, it is shown using scaling arguments that the decoupled fixed point is unstable 
against a fluctuation-induced biquadratic interaction, indicating that a crossover to nontrivial 
critical phenomena with different exponents appears as one approaches the critical point 
beyond the transient temperature region. This new scenario clearly contradicts the previous 
C ' notion of the 3D XY universality. 
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Introduction 



Field-induced critical phenomena of gapped spin-dimer antiferromagnets have drawn much 
attention. Such antiferromagnets typically consist of strongly coupled spin-1/2 dimers, and 
are essentially in singlet states in zero field. Elementary excitations in the gapped phase under 
external magnetic fields are S z = 1 triplet excitations, sometimes called "triplons," for which a 
magnetic field acts as a chemical potential. They undergo Bose-Einstein condensation (BEC) 
when their density is appropriately tuned. 1 ) 

BaCuSi20g is one of such spin-dimer compounds, 2 ~ 14 ) with characteristic frustration in 
interlayer interactions. Spin dimers in this compound align on the body-centered tetragonal 
(BCT) lattice (Fig. 1). Owing to the lattice geometry, loops that include interlayer hoppings 
are frustrated, which leads to cancellation among interlayer interactions. An interesting be- 
havior related to this interlayer frustration has been reported. 15 ^ The phase boundary around 
the quantum critical point is described by the power law T c (H) ex (H — H c )^ with an anoma- 
lous exponent <f> = l. 2 ) Since mean-field treatment yields <f> = 2/d, the exponent is regarded 
as a characteristic of two-dimensional (2D) systems, and in this sense the phenomenon is 
called "dimensional reduction." We refer to several recent papers for further details on this 
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Fig. 1. (Color online) (a) Localized spin degrees of freedom in BaCuSi206- (b) BCT lattice. For 
clarity, interlayer bonds are drawn partially. 

subject. 2 ' 10 - 14 ) 

Our main interest in this paper is to study the critical properties at the finite temperature 
transition in BaCuSi206- Owing to the broken U(l) symmetry of the ordered state and a few 
experimental observations such as the A-peak of specific heat, 6 ) it has been presumed so far 
that the transition is in the three-dimensional (3D) XY universality class. However, as we will 
argue in this paper, the additional Z2 symmetry breaking that characterizes the ordered phase 
in the BCT lattice makes the XY fixed point unstable. 10 ~ 13 ) The 3D XY-type order takes place 
in two subsystems or sublattices, namely, the even- and odd-numbered layers. The reason for 
this is that a Z2 mirror-image reflection symmetry of the BCT lattice precludes any bilinear 
effective coupling between layers on different sublattices. Although this Z2 symmetry allows 
biquadratic interlayer couplings making the order "collinear," whether the XY antiferromag- 
netic (AF) moments of the sublattices are parallel or antiparallel remains undetermined, and 
one of them is selected via spontaneous Z2 symmetry breaking. The structure of this symme- 
try breaking is clearly different from that for the standard XY ordering, and our final goal is 
to understand the effect of the additional Z2 symmetry breaking on the universality class of 
the transition. For this purpose, we will study a classical spin model that will be introduced 
in the next section. A Hamiltonian for classical spins is adequate for describing the critical 
behavior near the critical temperature because the relevant (largest) fluctuations are classical: 
fluctuations in the imaginary time direction become negligible because they are confined to 
a finite size (3 = 1/T. We will also present numerical simulation results that elucidate the 
above-mentioned symmetry-breaking structure. 

2. Model 

The low-energy subspace generated by the S z = 1 triplet and singlet well approximates the 
spin-dimer systems in an applied magnetic field. 1,16 ) The corresponding effective Hamiltonian 
is the XXZ model for 5 = 1/2 pseudo-spins that represent the two states of each dimer: 
S z = 1/2 for the triplet and S z = —1/2 for the singlet. Since the the thermodynamic phase 
transition is driven by classical (or thermal) fluctuations, we can replace S = 1/2 pseudo- 
spins with classical spins, in order to study the critical phenomena near the finite-temperature 
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transition: 

n,{r,r') n,r,Sr n,r 

where S ntr is a three-component classical spin, located at a site r in the n-th layer on the 
BCT lattice (r refers to a two-component vector), (r, r 1 ) are nearest-neighbor pairs on a given 
layer, and Sr = (±a/2, ±a/2) are interlayer displacement vectors (hereafter, we take a = 1). Jm 
( Jj_) is the AF intralayer (interlayer) interaction. Throughout the paper, we consider the case 
Jii > Jj_ > 0, whose inequality sign is the same as that in the relation between the magnitudes 
of interdimer exchanges in BaCuSi20g. The finite magnetic field H breaks 0(3) spin symmetry 
down to 0(2) symmetry, and thus the spins we treat are XY-like. This justifies the neglect of 
the easy-plane-type spin-anisotropy that exists in the effective XXZ model mentioned above. 
Although the two-component spins would serve the present purpose just as well as the three- 
component spins, we use the latter for technical reasons. 

First, we will discuss the basic properties of H c i. The energy is minimized when the spins 
form a canted AF order in each layer. The ground-state configuration is given by 

S n>r = (sin©cos$ n e iQ -( r - r o n) \ sin 9 sin ^ n e iQ < r ~ r o n) \ cos©), (2) 

with cos© = H/ [8(J)| + J±)\, Q = (vr,vr), and rj n) = (3, 5) <5(_i)» i- We define M xy = 
(cos sin $„) to represent the XY AF moment of the n-th layer, and in what follows we use 
the term "AF moment" to refer to this quantity unless otherwise specified. Interlayer mean- 
field interactions cancel out in the ground state because of a combination of the intralayer AF 
order and the geometry of the BCT lattice. This means that the AF moments of one layer 
can be rotated without changing ground-state energy. Thus, the system may be viewed as a 
set of independent 2D layers at T = 0. 

This apparent 2D character is lifted by thermal fluctuations. We will show that these 
fluctuations select a state qualitatively analogous to the ordered phase of the original quantum 
system. One of the simplest ways to see how this "order by disorder" 17,18 ) takes place is 
to use spin-wave approximation and evaluate free energy as a function of the ground-state 
configuration Let 9 n>r and 4>n,r represent small fluctuations around a given ground-state 

configuration. Expanding the Hamiltonian to the second order in these variables, we rewrite 
it in the form H c i « H sw (0, {$„} ; {#„,r}, {<f>n,r }) = "%2D + gV with g = J±/ J\\ . Here, 

= 9 ^ d ^ 0n,q0n,-q + ^ (<?) </>n,q<f>n,-q) (3) 
n,q 

and 

V = [^00 (9) 0n,q0n+l,-q + l^i ^ <f>n,q<j>n+l,-q + 7^ (q) (On,q<t>n+l,-q ~ 9n+l,-q<l> n ,q) , 
n,q 

(4) 

where 6 n , q = {L 2 y 1/2 £ r 6 n>r e~^ r and <j> n , q = {L 2 )~ 1/2 £ r <t>n,r ^ iq r with L 2 being the 
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number of sites in each layer. The coefficients in %2D are written as 

ujg (q) = 2 J|| [2 + (l — 2 cos 2 0) (cos q x + cos q y )] 
uj^ (q) = 2 J|| sin 2 6 (2 — cos q x — cos q y ) , 

and those in V are written as 
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(5) 
(6) 



1$ (q) = 4 J\\ [sin 2 6 cos f cos ^ - cos 2 9 cos ($ n+1 - $ n ) sin ^ sin ^ 



d (q) + C 2 (q) cos ($ n+ i - $ n ) 



7^ (q) 



-4 J|| sin 2 6 cos ($ n +i - $n) sin y sin ^ 
! l) (q) = 4 J|| sin G cos 9 sin ($ n+1 - $ n ) sin y sin ^, 



(7) 

(8) 
(9) 



where C\ = 4Jy sin 2 9 cos cos ^ and C2 = — 4Jy cos 2 9 sin ^ sin -y. Then we expand free 
energy as F = — TlnZo — TYX=i P k {y k ) < Q ) with 5 being a small parameter (Zo = 
Tre~^ 2D and (^ fc )g^ denote the cumulants with respect to Z ( ^~ 1 e _ ^^ 2D ). To the fourth order 
in g, we obtain 

2 4 
F = -T In Z - 9 -A (9) TL 2 £ cos 2 ($ n+1 - $ n ) - ^ S (9) TL 2 £ cos (<D n+2 -$„)-..., 



(10) 

with the definitions of A and B given below. We have dropped several "biquadratic" terms of 
0(g 4 ) because they do not change the 0{g 2 ) term's symmetry discussed below. 

The {$„}-dependence of F lifts part of the ground-state degeneracy, which is unrelated to 
the system symmetry. The coefficient A (9) of the 0(g 2 ) term is determined by the integral 



?<y*)f = i?Y,f 

„ Jbz 



d 2 q 



, . o "t" , ,0 +i" 



u e {q)- uj^qf co e (q)^(q) 
= A (9) L 2 ^2 cos2 ($n+i - $n) + const. 



(11) 



A is found to be positive, favoring collinear configurations where $ ra +i — <3? n = or ir. There- 
fore, this term represents the biquadratic effective interaction between nearest-neighbor lay- 
ers. On the other hand, there is no terms proportional to cos (& n +i — & n ) in F, being con- 
sistent with the fact that bilinear effective interactions between nearest-neighbor layers are 
forbidden by the symmetry of the BCT lattice. 10 ~ 13 ) A mirror-image transformation with 
respect to the (100) plane that contains (0,0) or (1/2,1/2) is such a symmetry operation. 
By the mirror-image transformation with respect to the plane that contains (0,0), for ex- 
ample, sites r with e 1 ^ '( r ~ r o > = ±1 in even- numbered (odd- numbered) layers are mapped 



,,, , , . „ v , . , -\ r - r o > = ±1 (e tv ' t - r-r o > = Tl)- Consequently, under 

the corresponding symmetry transformation S n ^ r — > S' nyr = S nt f, the local AF moments 
change as e^< r ~ r o n) ) S^ r -)• e i ^< r - r o n) ) S = (-l) n e iQ < f - r o n) )S^ (a = x,y). This means 
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$ n — > <J> n (n: even) and <!>„ — >■ + ir (n: odd) for the phase of their spatial average over 
a layer, resulting in cos ($ n +i — *^n) being mapped to — cos (<$n+i — 3>n)- Therefore, the two 
types of the collinear configurations, namely those with AF moments of adjacent layers being 
parallel (3> n +i — = 0) or antiparallel (|$ n+ i — <£ n | = 7r), are equivalent. While this de- 
generacy generally exists for the AF moments of any two layers (n,n') with \n — n'\ being an 
odd number, this is not the case with layers of even-numbered separations. Indeed, the 0(g A ) 
term in F determined by 



represents the bilinear effective interaction between second nearest-neighbor layers, favoring 
$«+2 - $„ = 0. 

Although the above spin-wave treatment describes the situation at low temperatures 
T <C T c , we can expect essentially the same form as eq. (10) also for T < T c in terms of 
symmetry. Therefore, the ordered phase is expected to have the following characteristics (see 
Fig. 2). First, there are two subsystems with XY-type 3D ordering, namely even- and odd- 
numbered layers, but there are no bilinear effective interactions in between. Second, because 
of the effective biquadratic interactions, the AF moments of these subsystems tend to align 
in the same direction. As a consequence, there are two symmetrically equivalent but distinct 
configurations 3>„+i — <£„ = or ir (Fig. 2). Note that interlayer bonds that are equivalent 
by symmetry become inequivalent in the ordered phase, meaning that bond ordering results 
from the spontaneous Z2 symmetry breaking. These features are qualitatively the same as the 
original quantum system. 10 ~ 13 ) The bond order is a direct 3D analogue of the "Ising-order" 
that was discussed in the frustrated square-lattice J1-J2 model for U2 > Ji. 18-20 ) We will use 



= E { S n,rS X n+ l,r + Sr + ^M+l.r+ft.) - ("!)* = ^ [<(*, -*) " H (13) 



as the local bond-ordering order parameter. It is invariant under 0(2) spin rotations but 
changes its sign (a — > — a) under mirror-image reflections of the lattice, i.e., it serves to detect 
the Z2 symmetry breaking. 

3. Theoretical Arguments on the Phase Transition 

3.1 Single bilayer and the BCT lattice 

Let us consider the problem of how many phase transitions take place. There are two 
possible scenarios: a) only one phase transition at T = T^ Y = T c BO , driven by the XY 
ordering, b) two transitions at T = Tf° and T = T* Y (< T c BO ). Since the XY ordering 
necessarily accompanies the bond ordering, T* Y > T®° is impossible. We will use the results 
of a single bilayer case to argue that a) is the correct scenario. The bilayer case is equivalent to 
the J1-J2 XY model with J\ = J± and J2 = J\\. This model has been studied numerically by 
Loison and Simon. 20 ) It breaks the additional Z2 symmetry for 2 J2 > Ji at a finite temperature 




(12) 
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Fig. 2. Expected ordered phase. Solid and gray arrows denote spins on different sublattices, namely, 
even- and odd-numbered layers, and only their XY components are shown for clarity. These 
configurations are distinct in the sense of the spontaneous bond ordering (see text). Solid (dashed) 
interlayer lines represent the spin-pairs with parallel (antiparallel) XY components. 



T = T C K with a second-order transition, which is followed by a Berezinskii-Kosterlitz-Thouless 
(BKT)-type transition at a slightly lower temperature T = Tkt < T£. 

We next consider an array of weakly coupled bilayers with Jf being the inter-bilayer 
interaction (Jy/Jn <C 1), such that the array returns to the original BCT lattice when Jf = J±. 
In this case, induced by the order-by-disorder mechanism, there appear effective interlayer 
interactions J' R g = (Jy / J\\) k J\\ between the XY-components of spins and effective interlayer 
interactions J"ff= {Jf/J\\f J\\ between bond-ordering order parameters. Equation (10) implies 
that k = 4, and also that I = 2 because J"g is determined by effective biquadratic interactions. 
The effective coupling J' e g is relevant for turning the BKT transition into the true long-range 
XY ordering and J" e ^ drives the 2D bond ordering to the 3D behavior. To determine T XY and 
T c BO , we use a simple random phase approximation (RPA) argument. 21 ) In this treatment, 
the XY and bond orderings take place when 



Jt 



/r\ I',:" J ~ i (14) 
J'I ff X K (T c BO ) « 1 (15) 

are satisfied, respectively. Here, x (T) and x K (T) are the AF XY and 2D bond ordering suscep- 
tibilities for the single bilayer, respectively. As the temperature is lowered, x (T) is expected 
to diverge exponentially as x(T) oc J|| _1 exp (byj jr^^pj with b being a constant. On the 

other hand, x K (T) is expected to show the power-law divergence x K (T) oc Jy ( TK C J 
(7 > 0). By substituting these expressions in eqs. (14) and (15), we obtain: 

2 

(16) 



T C XY - T KT ^ 


b 


2 


b/k 


Tkt 


Aw). 




_ln(j//J||)_ 



rpBO 

c 



rpK 
C 



{J%U\\ 



,1/7 



(Jf/J\\ 



(17) 



Because eq. (16) diverges as Jf approaches Jn while eq. (17) does not, these equations imply 
that T c xy > T c BO for J f > Jf. Here, Jf depends on the difference T£-T KT . Since the difference 
seems to be very small according to the existing numerical simulations, 22 ) we can expect that 
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J'j is small. Therefore, we can conclude that the thermodynamic phase transition of the XY 
ordering first takes place as temperature decreases over a wide range of J/ values. Here, note 
that the above RPA estimate of is based on assumption that the XY spin ordering 

is absent. Since the XY ordering also breaks the Z2 symmetry, the bond ordering transition 
temperature cannot be lower than T^ Y . Therefore, the above RPA result T c xy > Tf° actually 
implies a single-phase transition. 

3.2 Stability of the decoupled XY fixed point 

The next question is about the universality class of the phase transition, in particular as 
to whether it belongs to the previously expected 3D XY universality class. Introducing the 
"continuous" 0(2) real vectors (f>f(r) (a = x,y) to describe the spins on the even- (i = 1) 
and odd-numbered (i = 2) layers, 23 ) we consider the Landau-Ginzburg- Wilson (LGW)-type 
effective Hamiltonian of the form 

H *8 = f [I E ■ Mi + ^2 • d^ 2 ) + t (|^| 2 + |<H 2 ) + u (|0i| 4 + \<P 2 \ A ) 

+ A(0i ■ cp 2 f + g\<Pi\ 2 \h\ 2 } d d r. (18) 

The first three terms constitute a standard cp 4 theory for the decoupled 0(2) model. The 
(<Ai • 02) 2 term represents a quadrupole-quadrupole interaction induced by the order-by- 
disorder mechanism. The other quartic term |^i| 2 |^2| 2 is included here explicitly, because 
it is generated through renormalization group (RG) iterations. As we mentioned earlier, the 
lattice-symmetry of the original model eq. (1) does not allow effective bilinear interactions 
between nearest-neighbor layers. For this reason, the quadratic term 4>\ ■ <p 2 is impossible in 
eq. (18). 

This Hamiltonian eq. (18) is an iV = M = 2 case of the model referred to as the "iV-coupled 
M-vector model," 24 ) with the additional ((pi ■ fa) 2 term. It has a trivial decoupled fixed point 
(D) at u 7^ and A = g = 0, which is a plausible candidate for the fixed point corresponding 
to the expected 3D XY universality. This model was first introduced in the 1970s, 25 ) in the 
context of the replica theory for random systems. It was found that the decoupled fixed point 
is unstable against perturbations such as biquadratic ones. Below, we briefly summarize the 
argument, because the original argument was made in a relatively different context. 

The stability of a fixed point against a given perturbation is determined by the RG of its 
conjugate field. Therefore, we need to compute the RG eigenvalues y\ t D {Ug,D) of the coupling 
A (g) at the decoupled fixed point. They can be computed via two-point correlators at the 
decoupled fixed point, which in this case can be readily factorized into known correlators. 
First, 

<|<Mr)| 2 |<Mr)| 2 IMrOl'l'MrOl^ = {\^{r)\ 2 \Mr')?) D <l<Mr)| 2 |<Mr')| 2 ) D 
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cx |r - r'r 4x *, (19) 

where xt = d — 1/u is the scaling dimension of the energy-density field of the 3D XY model 
with v being the correlation- length exponent. This means that the scaling dimension of the 
101 1 2 102 1 2 term is equal to 2x t . Therefore, 

VgD = d-2x t = 2/v-dxi -0.021815 < in d = 3, (20) 

where we used v = 0.67155(27). 26 ) The negative y 9t D indicates that the decoupled fixed point 
is stable against the |0i| 2 |</>2| 2 term. However, this is not the case with the other (cfti ■ 4>2) 2 
term. This term can be decomposed into 

(<k • 2 ) 2 = \ (QTQT + Q7Q7 + \M 2 \M 2 ) , (21) 

where Qf x = (<fif) 2 — (4>f) 2 and Q* y = 2(f>f(f> y (i = 1,2) are components of the traceless 
symmetric tensor of the quadrupole order parameter. Using eq. (21) and the 0(2) invariance 
of H e ff, we obtain 

4<(0i-0 2 ) 2 (r) ^i-fo) 2 {r')) D 

= 2 {QT{r)QT(r')) D {QT(r)Q'T(r')) D + <0 2 (r)0 2 (r')} D (^(r)0 2 (r')> D 

^ Cqq C u (22 , 

| r _ r l\AxQ | r _ r l\Axt ' ^ ' 

where Cqq and Cu are nonzero coefficients and xq is the scaling dimension of the quadrupole 
order parameter. Comparing xq 1.237 27 ) with x t ~ 1.5 1 09, 26 ) we find that the quadrupole- 
quadruple correlator gives the most relevant contribution to eq. (22). Consequently, the scaling 
dimension of the (4>\ ■ fo) 2 term is equal to 2xq and we obtain 

y x ,D = d-2x Q ^ 0.526. (23) 

The positive y\ t D indicates that A is a relevant coupling for the decoupled 3D XY fixed point. 
In other words, the decoupled fixed point is unstable under such perturbation. 

4. Results of the Monte Carlo Simulation 



In this section, we present the results our Monte Carlo (MC) simulation for the Hamil- 
tonian eq. (1). All the results shown in what follows are obtained for J±/J\\ = 0.75 and 
H/J\\ = 10.0. The system size is L x L x (L/4) with L = 16,24,32,40, and 48. Our main 
motivation for these parameters and the anisotropic aspect ratio is to realize the proper 
configuration in finite-size systems. These considerations are necessary, because fluctuation- 
induced effective couplings, if normalized per site, are typically smaller than Jy by 2-3 orders 
of magnitude. 

Figures 3 and 4 show low-temperature snapshots of MC simulations, which are useful for 
understanding the situation. The points in these figures represent the local XY AF moments 
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Fig. 3. Low-temperature snapshot of a MC simulation for L = 48 at T/J\\ = 0.27 (see text). 
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Fig. 4. Snapshot under the same conditions as the snapshot in Fig. 3, but obtained using a different 
random-number sequence. 



e iQ-{r-r[ ) ){S%, Sn) in each layer, with their spatial positions being discarded for clarity; the ar- 
rows represent the in-layer AF moments. We can see that the configurations are collinear, and 
that the AF-moments in every other layer tend to align in the same direction. By comparing 
these two figures, we can also see that the effective interaction between layers of odd-number 
separations is not bilinear but biquadratic. As we will quantitatively show below, these con- 
figurations are typical at low temperatures. 

In Figs. 5(a) and 5(b), we show decays of the correlation function G(R) = 
G (R(n,r;n',r')) = (S% r S*, r , + Sn, r S^, r ,) along the interlayer [111] and intralayer [110] 
directions. G(R) along [110] indicates the formation of the in-layer AF order, and positive 
correlations along [111] for R/(l/2, 1/2, 1/2) = 2 and 4 suggest the presence of a bilinear 
effective interlayer interaction of the ferromagnetic type. On the other hand, the suppression 
of G(R) for R/ (1/2, 1/2, 1/2) = 1 and 3 is mainly due to the cancellation of the contributions 
of opposite ("parallel" and "antiparallel" ) configurations, suggesting the absence of a bilinear 
effective interaction between nearest-neighbor layers; nonzero values at these distances are 
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Fig. 5. Correlation function G(R) along [111] (a) and [110] (b) for L = 32. The lines are guides to 
the eyes. 



due to short-wavelength fluctuations. 

Next we turn to the analysis of the order parameters. There are two relevant order pa- 
rameters to be examined: the XY and bond-ordering order parameters. We use 

Mxy = £ e^(—o n) ) (Sl r x + S» r y) , (24) 

ngcvcn, r 

as the XY order parameter, namely, an AF order parameter defined on the "even" sub- 
lattice (the choice of either "even" or "odd" is arbitrary). On the other hand, we use 
Mbo = ^2 nr a n,r as the bond-ordering order parameter (fT n ,r is defined by eq. (13)). 

Figures 6(a) and 6(b) show the temperature dependence of the Binder parameters U^ Y = 
(M XY )/(M XY ) 2 and Uf° = (M BO )/(M 2 BO ) 2 . They should asymptotically show crossings 
for different system sizes at critical points, whereas U± Y — > 2 and Uf° — > 3 for T/Jy 3> 1, 
and U? Y , Uf° 1 for T / ' J\\ < 1. It is not easy to determine T^ Y and precisely as they 
still suffer from severe finite-size effects, but it is clear that there is a phase transition. Both 
Tf Y and Tf° are located in the region 0.305 < T/Jy < 0.310. 

In order to obtain the critical exponents, we perform finite-size scaling analysis of the 
squared quantities (M XY ) and {^bo) assuming the following standard scaling forms: 

(M 2 XY ) = L- (d - 2+r » ^ XY (tL^), (25) 

(M 2 BO ) = L~ 2 ^ ^Boit'L 1 ^'), (26) 

where t = (T - T c xy ) /T c xy and t' = (T - Tf°) )Tf° are reduced temperatures, $xy and 
&bo are scaling functions, x c is the scaling dimension of the bond order parameter, and the 
exponents v, v', and rj are conventional parameters. As shown in Figs. 7(a) and 7(b), we can 
produce a reasonable data collapse, where we use 

7] = 0.04(1), v = 0.67(1), and T c xy /J,| = 0.308(1) (27) 
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Fig. 6. Temperature dependence of the Binder parameters of the (a) XY ordering and (b) bond 
ordering. The insets show enlarged views in the critical region. The lines are guides to the eyes. 



for the XY ordering and 

2x a = 2.07(1), v' = 0.67(1), and T C B ° / J\\ = 0.310(2) (28) 

for the bond ordering. These values of T^ Y and are consistent with the estimations 

made from the Binder parameters. 

Although there is a slight difference between T^ Y and , the close proximity of v and 
v' (no difference in our resolution) suggests that the difference found for L < 48 is due to 
finite-size effects and that there is only one transition. As for the critical exponents, their 
values agree with those of the decoupled 3D XY fixed point, for which t]d and coincide 
with those of the 3D XY model (77 = 0.0380(4) and v = 0.67155(27) 26 ) ) and, in addition, 
x c ,D = 2x holds with x being the scaling dimension of the XY order parameter (therefore, 
2x^ D = 2 • 2x = 2(d - 2 + 77) = 2.0760(8) ). 13) The reason for this is that a corresponds to 
(f>i ■ <p2 in the continuous-spin language. Its two-point correlator can be factorized as 

<(</>! • (r) (0i • <h) {r')) D = (Mr) ■ Mr')) D {<h{r) ■ Mr')) D (29) 

using the 0(2) symmetry of the Hamiltonian eq. (1) or (18). Then, the simple counting of 
powers leads to x a> D = 2x. 
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Fig. 7. Finite-size scaling plots of (a) (M^y) and (b) (M| ). 



5. Discussion 

The MC results apparently contradict the previous observation that the decoupled fixed 
point is unstable, but the reason why critical exponents of such an unstable fixed point are 
obtained can be understood through the theory of crossover behavior. First, the above MC 
results suggest that the actual RG flow passes through the vicinity of the decoupled fixed 
point. In this parameter region, the system-size dependence of the singular part of free energy 
has the scaling form 

f s (t,\,L) = \t\ dUD ^> (tL 1 /" D , AL^.°) . (30) 

The amplitude of the relevant biquadratic coupling A can be estimated as |A| ~ ( J±/J\\) 2 A 
0.022 using eq. (10). It is too small to observe any significant sign of the crossover behavior in 
the sense that we need a system as large as L ~ |A| _1 / S/A - D 1400, which is much larger than 
the largest system investigated in the present work. In this way, the most natural interpretation 
of our MC results is summarized as follows: what we obtained are "effective" exponents due 
to transient behavior of the RG flow around the decoupled fixed point. However, since the 
flow should eventually leave the vicinity of the decoupled fixed point, we expect a crossover 
to be observed in larger systems. 

The above finite-size scaling arguments can be immediately translated into general scaling 
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arguments, in which experimental relevance becomes more transparent. The fact that the RG 
flow passes through the vicinity of the decoupled fixed point means that there is a certain 
temperature region where long-wavelength fluctuations are described by the decoupled 3D 
XY model. However, because this fixed point is unstable, this temperature region must be 
outside of the ultimate critical region: \t\ > tx such that the correlation length is bounded as 
f < lAI-VfA.D.as) We expect that, as one approaches the critical point beyond this intermediate 
temperature region, a crossover from the decoupled 3D XY behavior appears. 

A question that now arises concerns the nature of the stable fixed point. Since the pertur- 
bation field A is only slightly relevant at the decoupled fixed point for d sufficiently close to 4 
(y\,D = | € + 0(e 2 ) with e = 4 — d), it does not seem to be difficult to find a new stable fixed 
point by performing the e-expansion in the lowest nontrivial order. However, the attempt to 
conduct 0(e) calculation is unsatisfactory. 25 ) According to such a calculation, there is no sta- 
ble fixed point in the finite region of the parameter space. Moreover, it fails to reproduce the 
irrelevance of g around the decoupled fixed point, implying that the topology of the flow itself 
may not be convincing unless one goes to sufficiently high orders in e. Therefore, it remains to 
be clarified whether a stable fixed point exists and, if it does, what type of critical behavior is 
expected from such a fixed point. We expect that direct numerical simulations of the effective 
Hamiltonian eq. (18) will shed light on this problem. Such a simulation will be performed in 
a future study. 

6. Conclusions 

In this paper, we have studied a classical model of the finite-temperature transition in 
BaCuSi2C>6. We have demonstrated that thermal fluctuations select a particular configura- 
tion in a low-temperature ordered phase. It is stabilized by a composite of two subsystems, 
and is realized via multiple symmetry breakings, 0(2) and Z2, with the latter being related 
to the mirror-image reflection symmetry of the underlying lattice. The qualitative character- 
istics of the phase are the same as those of the original quantum system. On the basis of 
the RPA argument, we have also argued that there is only one transition at which the XY 
and bond orderings occur simultaneously in contrast to the 2D single-bilayer case. As for 
the critical behavior of the phase transition, we have shown that a plausible form of LGW- 
type effective Hamiltonian has a decoupled fixed point that yields the critical exponents of 
the 3D XY universality class. However, scaling arguments reveal that the decoupled fixed 
point is unstable against the perturbation of the quadrupole-quadrupole interaction. Since 
the quadrupole-quadrupole interaction is induced by the order-by-disorder effect, this con- 
clusion clearly excludes the decoupled XY fixed point for describing the asymptotic critical 
behavior. On the other hand, MC simulation yields a set of exponents that can be attributed 
to the decoupled 3D XY universality. As a first point to be noticed from this observation, we 
have indicated that there is actually an intermediate temperature region outside the ultimate 
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critical region, where the observed thermodynamic properties are related to the decoupled 
fixed point. We have also presented a reasonable argument explaining why no significant devi- 
ation has been observed (within our system sizes) from the critical exponents of the unstable 
decoupled fixed point, on the basis of the theory of crossover behavior. 
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